Project 9
CodesRUS
12/13/2021
Introduction
Trans-gulf migrants, known as TGMs, are species that endure the long migration over the Gulf of Mexico from Central and South America before arriving on land between Texas and Florida. The birds usually migrate during the springtime in a nonstop flight in hopes to find newer land in which they can forage and breed on. The specific time that the birds arrive on the North American land can have major impacts on the fitness of the bird overall. Birds who arrive sooner tend to have better opportunities to find a mate, whereas birds that arrive later have a more ensured food stability. In response to climate change and warming temperatures, TGMs need to shift the date in which they arrive at their breeding grounds or we could see mass population declines. Scientists are currently debating the extent to which TGMs and other long-distance migrants are shifting their arrival patterns, if they are at all. In order to determine the magnitude of these shifts, scientists have been using field, tracking, and banding methods but data has been limited as there are not many scientists who do this type of work.
In more recent years, trends over species conservation and climate change awareness have inspired further research and developments in similar issues. One such study conducted by (Charmantier et al. (2008)) highlighted the importance of understanding the mechanisms behind population responses to climate change. Phenotypic plasticity is one such mechanism they studied and found that individual responses to these environmental changes may be fixed in the population but evolve over time. Specifically concerning TGMs, many online databases have arisen that allow you to log observations about birds that can be seen by the public. One in particular, fittingly titled eBird, is one of the most complex and has proven to be especially helpful to scientists. Created by the Cornell Laboratory of Ornithology and the National Audubon Society, eBird launched in 2002 and has over 100 million bird observations made by thousands of different people. The database even allows observers to log historical data so observations from pre-2002 can be read as well. In this project we want to use eBird and meteorological data to study what affects the local weather conditions can have on the arrival times of TGMs in the state of Massachusetts.
Methods
In order to find the effect of local weather conditions on the arrival time of TGMs in Massachusetts, we first used a for loop to load all the data of human observation in Massachusetts from eBird database. Then, we used a token obtained from NOAA website to download weather data of three sample locations. After that, we moved to the data analysis section. To prepare the eBird data, we plotted the proportion of the population that arrived and predicted their arrival time. At the end of this part, we plotted how this arrival day varies with year to see the influence of climate on arrival date. To prepare for weather data, we calculated wind direction and joined the weather data with the species data. Next, we calculated the mean of our weather variables 5 and 10 days beforehand by using the frollmean() function to get more accurate information. Finally, we used lme from lme4 package to produce Linear Mixed-effect Model. Anova tested were performed to find the best fit model and dredg() function was used for model testing.
Results
### choice of species
#Chimney Swift, Chaetura pelagica
#Ruby-throated Hummingbird, Archilochus colubris
#Belted Kingfisher, Megaceryle alcyon
#Yellow-bellied Sapsucker, Sphyrapicus varius
#Common Yellowthroat,Geothlypis trichas
swift <- occ_data(scientificName = "Chaetura pelagica", stateProvince = "Massachusetts", limit = 200, year = 2018)
humming<- occ_data(scientificName = "Archilochus colubris", stateProvince = "Massachusetts", limit = 200, year = 2018)
kingfisher<- occ_data(scientificName = "Megaceryle alcyon", stateProvince = "Massachusetts", limit = 200, year = 2018)
sapsucker<- occ_data(scientificName = "Sphyrapicus varius", stateProvince = "Massachusetts", limit = 200, year = 2018)
yellowthroat<- occ_data(scientificName = "Geothlypis trichas", stateProvince = "Massachusetts", limit = 200, year = 2018)
MA <- map_data('state', 'massachusetts')
swift.p <- ggplot(MA, aes(long,lat,group=subregion) )+
geom_polygon(colour = "gray",fill="gray90")+geom_point(data=swift[[2]],aes(x=decimalLongitude,y=decimalLatitude,size=individualCount),alpha=0.3,inherit.aes = F)+ coord_quickmap()+theme_void()
humming.p <- ggplot(MA, aes(long,lat,group=subregion) )+
geom_polygon(colour = "gray",fill="gray90")+geom_point(data=humming[[2]],aes(x=decimalLongitude,y=decimalLatitude,size=individualCount),alpha=0.3,inherit.aes = F)+ coord_quickmap()+theme_void()
kingfisher.p <- ggplot(MA, aes(long,lat,group=subregion) )+
geom_polygon(colour = "gray",fill="gray90")+geom_point(data=kingfisher[[2]],aes(x=decimalLongitude,y=decimalLatitude,size=individualCount),alpha=0.3,inherit.aes = F)+ coord_quickmap()+theme_void()
sapsucker.p <- ggplot(MA, aes(long,lat,group=subregion) )+
geom_polygon(colour = "gray",fill="gray90")+geom_point(data=sapsucker[[2]],aes(x=decimalLongitude,y=decimalLatitude,size=individualCount),alpha=0.3,inherit.aes = F)+ coord_quickmap()+theme_void()
yellowthroat.p <- ggplot(MA, aes(long,lat,group=subregion) )+
geom_polygon(colour = "gray",fill="gray90")+geom_point(data=yellowthroat[[2]],aes(x=decimalLongitude,y=decimalLatitude,size=individualCount),alpha=0.3,inherit.aes = F)+ coord_quickmap()+theme_void()
swift.p2 <- ggdraw() + draw_image("swift.png", scale =0.3, halign=0, valign=1) + draw_plot(swift.p)
print(swift.p2)
humming.p2 <- ggdraw() + draw_image("humming.png", scale =0.3, halign=0, valign=1) + draw_plot(humming.p)
print(humming.p2)
kingfisher.p2 <- ggdraw() + draw_image("kingfisher.png", scale =0.3, halign=0, valign=1) + draw_plot(kingfisher.p)
print(kingfisher.p2)
sapsucker.p2 <- ggdraw() + draw_image("sapsucker.png", scale =0.3, halign=0, valign=1) + draw_plot(sapsucker.p)
print(sapsucker.p2)
yellowthroat.p2 <- ggdraw() + draw_image("yellowthroat.png", scale =0.3, halign=0, valign=1) + draw_plot(yellowthroat.p)
print(yellowthroat.p2)
The above are the records for each of the 5 species selected in Massachusetts, 2018.
species <- c("Chaetura pelagica","Archilochus colubris","Megaceryle alcyon","Sphyrapicus varius","Geothlypis trichas")
y <- paste0("2000",",","2019")
m <- paste0("4",",","5")
dat.l <-list()
for(s in species){
n.obs <- occ_data(scientificName = s,year=y,month=m,limit=0,country="US",basisOfRecord = "HUMAN_OBSERVATION",stateProvince="Massachusetts")$meta$count
print(n.obs)
dat.l[[paste0(s)]] <- occ_data(scientificName = s,year=y,month=m,
limit=n.obs,country="US",
basisOfRecord = "HUMAN_OBSERVATION",
stateProvince="Massachusetts")[[2]]
}
dat <- rbindlist(dat.l,fill=T)
head(dat)
saveRDS(dat,"massbird.data.RDS")
library(tidyverse)
dat <- readRDS("massbird.data.RDS")
dat.l%>%
bind_rows(.id = "species") %>%
group_by(year,species)%>%
summarise(count=sum(individualCount,na.rm = T))%>%
ggplot(aes(x=year,y=count,col=species))+geom_point()
To begin our analysis, collection and sorting of data needed to occur for each species selected within a specific state during a range of years. A data set was created to do such.
options(noaakey = "VYpNYgvfxWWfUIWmcjVZxUUJwebzTmLB")
sts <- c(
"GHCND:USW00013894", #Mobile, AL 2k away about 10 days away @200 km/day
"GHCND:USW00013881", #Charlotte, NC 1000 km away about 6 days away @200 km/day
"GHCND:USW00014739" #Boston
)
bos <- ncdc_stations(stationid = "GHCND:USW00014739")
print(bos)
cha <- ncdc_stations(stationid = "GHCND:USW00013881")
print(cha)
mob <- ncdc_stations(stationid = "GHCND:USW00013894")
print(mob)
sta.d <- bind_rows( #bind the rows
lapply(sts,function(x) ncdc_stations(stationid = x)$data )
)%>%
left_join(usmap_transform(.[,c("longitude","latitude")]))%>%
mutate(name=str_sub(name, -5,-4))%>%
mutate(migr.day=c(10,5,0))%>%
separate(id,into = c("station.type","id"))%>%
print()## Joining, by = c("latitude", "longitude")plot_usmap(
include = c(.northeast_region,.south_region,.east_north_central)
)+geom_point(data=sta.d,aes(x=longitude.1,y=latitude.1,col=name),size=5)+geom_label(data=sta.d,aes(x=longitude.1,y=latitude.1,col=name,label=name),size=5,nudge_x = 1e6*0.25)+theme(legend.position = "none")
weather.d <- meteo_pull_monitors(sta.d$id,date_min = "2005-01-01") #since 2005 we see an uptick in eBird data## using cached file: ~/Library/Caches/R/noaa_ghcnd/USW00013894.dly## date created (size, mb): 2021-12-22 17:31:06 (7.048)## file min/max dates: 1948-01-01 / 2021-12-31## using cached file: ~/Library/Caches/R/noaa_ghcnd/USW00013881.dly## date created (size, mb): 2021-12-22 17:31:24 (7.645)## file min/max dates: 1939-07-01 / 2021-12-31## using cached file: ~/Library/Caches/R/noaa_ghcnd/USW00014739.dly## date created (size, mb): 2021-12-22 17:31:45 (8.477)## file min/max dates: 1936-01-01 / 2021-12-31head(weather.d)Using NOAA’s NCDC API, weather parameters during certain time frames were collected and placed within data sets. Weather sampling locations for this project are represented above.
cs<- dat.l%>%
bind_rows(.id = "species") %>%
filter(species=="Chaetura pelagica")%>%
group_by(year)%>%
mutate(date=as.Date(paste0(year,"-",month,"-",day)),
j.day=julian(date,origin=as.Date(paste0(unique(year),"-01-01")))
)%>%
group_by(species,year,j.day,date)%>%
summarise(day.tot=sum(individualCount,na.rm=T))%>%
group_by(species,year)%>%
mutate(prop=cumsum(day.tot/sum(day.tot,na.rm = T)))%>%
filter(year>2004)## `summarise()` has grouped output by 'species', 'year', 'j.day'. You can override using the `.groups` argument.cs%>%
ggplot(aes(j.day,prop))+geom_point()+facet_wrap(year~.)
cs.pred <- cs%>%
group_by(year)%>%
summarize(
pred=predict(nls(prop~SSlogis(j.day,Asym, xmid, scal)),newdata=data.frame(j.day=min(j.day):max(j.day))),
j.day=min(j.day):max(j.day),
)%>%
left_join(cs%>%select(j.day,date)) ## `summarise()` has grouped output by 'year'. You can override using the `.groups` argument.## Adding missing grouping variables: `species`, `year`## Joining, by = c("year", "j.day")cs%>%
ggplot(aes(j.day,prop))+geom_point(aes=0.3)+geom_line(data=cs.pred,aes(x=j.day,y=pred),col="blue",size=2)+facet_wrap(year~.)
cs.arrive.date <-cs.pred%>%
group_by(year)%>%
filter(j.day==j.day[which.min(abs(pred-0.25))])
cs.arrive.date%>%
ggplot(aes(year,j.day))+geom_point()
For Chimney Swift, the arrivals are depicted above with the estimation of arrivals. Then, using this data we plotted how the arrival day varied depending on the year. Based on the plots, the Julian day when 25% of the population arrived were 127-130. For all species, there was no clear cut trend that occurred in arrival date as time continued.
rh<- dat.l%>%
bind_rows(.id = "species") %>%
filter(species=="Archilochus colubris")%>%
group_by(year)%>%
mutate(date=as.Date(paste0(year,"-",month,"-",day)),
j.day=julian(date,origin=as.Date(paste0(unique(year),"-01-01")))
)%>%
group_by(species,year,j.day,date)%>%
summarise(day.tot=sum(individualCount,na.rm=T))%>%
group_by(species,year)%>%
mutate(prop=cumsum(day.tot/sum(day.tot,na.rm = T)))%>%
filter(year>2004)## `summarise()` has grouped output by 'species', 'year', 'j.day'. You can override using the `.groups` argument.rh%>%
ggplot(aes(j.day,prop))+geom_point()+facet_wrap(year~.)
rh.pred <- rh%>%
group_by(year)%>%
summarize(
pred=predict(nls(prop~SSlogis(j.day,Asym, xmid, scal)),newdata=data.frame(j.day=min(j.day):max(j.day))),
j.day=min(j.day):max(j.day),
)%>%
left_join(rh%>%select(j.day,date)) ## `summarise()` has grouped output by 'year'. You can override using the `.groups` argument.## Adding missing grouping variables: `species`, `year`## Joining, by = c("year", "j.day")rh%>%
ggplot(aes(j.day,prop))+geom_point(aes=0.3)+geom_line(data=rh.pred,aes(x=j.day,y=pred),col="blue",size=2)+facet_wrap(year~.)
rh.arrive.date <-rh.pred%>%
group_by(year)%>%
filter(j.day==j.day[which.min(abs(pred-0.25))])
rh.arrive.date%>%
ggplot(aes(year,j.day))+geom_point()
For Ruby-throated Hummingbird, the Julian day when 25% of the population arrived were 129-132.
bk<- dat.l%>%
bind_rows(.id = "species") %>%
filter(species=="Megaceryle alcyon")%>%
group_by(year)%>%
mutate(date=as.Date(paste0(year,"-",month,"-",day)),
j.day=julian(date,origin=as.Date(paste0(unique(year),"-01-01")))
)%>%
group_by(species,year,j.day,date)%>%
summarise(day.tot=sum(individualCount,na.rm=T))%>%
group_by(species,year)%>%
mutate(prop=cumsum(day.tot/sum(day.tot,na.rm = T)))%>%
filter(year>2004)## `summarise()` has grouped output by 'species', 'year', 'j.day'. You can override using the `.groups` argument.bk%>%
ggplot(aes(j.day,prop))+geom_point()+facet_wrap(year~.)
bk.pred <- bk%>%
group_by(year)%>%
summarize(
pred=predict(nls(prop~SSlogis(j.day,Asym, xmid, scal)),newdata=data.frame(j.day=min(j.day):max(j.day))),
j.day=min(j.day):max(j.day),
)%>%
left_join(bk%>%select(j.day,date)) ## `summarise()` has grouped output by 'year'. You can override using the `.groups` argument.## Adding missing grouping variables: `species`, `year`## Joining, by = c("year", "j.day")bk%>%
ggplot(aes(j.day,prop))+geom_point(aes=0.3)+geom_line(data=bk.pred,aes(x=j.day,y=pred),col="blue",size=2)+facet_wrap(year~.)
bk.arrive.date <-bk.pred%>%
group_by(year)%>%
filter(j.day==j.day[which.min(abs(pred-0.25))])
bk.arrive.date%>%
ggplot(aes(year,j.day))+geom_point()
For Belted Kingfisher, the Julian day when 25% of the population arrived were 105-108.
ys<- dat.l%>%
bind_rows(.id = "species") %>%
filter(species=="Sphyrapicus varius")%>%
group_by(year)%>%
mutate(date=as.Date(paste0(year,"-",month,"-",day)),
j.day=julian(date,origin=as.Date(paste0(unique(year),"-01-01")))
)%>%
group_by(species,year,j.day,date)%>%
summarise(day.tot=sum(individualCount,na.rm=T))%>%
group_by(species,year)%>%
mutate(prop=cumsum(day.tot/sum(day.tot,na.rm = T)))%>%
filter(year>2004)## `summarise()` has grouped output by 'species', 'year', 'j.day'. You can override using the `.groups` argument.ys%>%
ggplot(aes(j.day,prop))+geom_point()+facet_wrap(year~.)
ys.pred <- ys%>%
group_by(year)%>%
summarize(
pred=predict(nls(prop~SSlogis(j.day,Asym, xmid, scal)),newdata=data.frame(j.day=min(j.day):max(j.day))),
j.day=min(j.day):max(j.day),
)%>%
left_join(ys%>%select(j.day,date)) ## `summarise()` has grouped output by 'year'. You can override using the `.groups` argument.## Adding missing grouping variables: `species`, `year`## Joining, by = c("year", "j.day")ys%>%
ggplot(aes(j.day,prop))+geom_point(aes=0.3)+geom_line(data=ys.pred,aes(x=j.day,y=pred),col="blue",size=2)+facet_wrap(year~.)
ys.arrive.date <-ys.pred%>%
group_by(year)%>%
filter(j.day==j.day[which.min(abs(pred-0.25))])
ys.arrive.date%>%
ggplot(aes(year,j.day))+geom_point()
For Yellow-bellied Sapsucker, the Julian day when 25% of the population arrived were 105-113.
#Common Yellowthroat,Geothlypis trichas(referred by cn)
cn<- dat.l%>%
bind_rows(.id = "species") %>%
filter(species=="Geothlypis trichas")%>%
group_by(year)%>%
mutate(date=as.Date(paste0(year,"-",month,"-",day)),
j.day=julian(date,origin=as.Date(paste0(unique(year),"-01-01")))
)%>%
group_by(species,year,j.day,date)%>%
summarise(day.tot=sum(individualCount,na.rm=T))%>%
group_by(species,year)%>%
mutate(prop=cumsum(day.tot/sum(day.tot,na.rm = T)))%>%
filter(year>2004)## `summarise()` has grouped output by 'species', 'year', 'j.day'. You can override using the `.groups` argument.cn%>%
ggplot(aes(j.day,prop))+geom_point()+facet_wrap(year~.)
cn.pred <- cn%>%
group_by(year)%>%
summarize(
pred=predict(nls(prop~SSlogis(j.day,Asym, xmid, scal)),newdata=data.frame(j.day=min(j.day):max(j.day))),#predict the logistic curve for each species
j.day=min(j.day):max(j.day),
)%>%
left_join(cn%>%select(j.day,date))## `summarise()` has grouped output by 'year'. You can override using the `.groups` argument.## Adding missing grouping variables: `species`, `year`## Joining, by = c("year", "j.day")cn%>%
ggplot(aes(j.day,prop))+geom_point(aes=0.3)+geom_line(data=cn.pred,aes(x=j.day,y=pred),col="blue",size=2)+facet_wrap(year~.)
cn.arrive.date <-cn.pred%>%
group_by(year)%>%
filter(j.day==j.day[which.min(abs(pred-0.25))])
cn.arrive.date%>%
ggplot(aes(year,j.day))+geom_point()
For Common Yellowthroat, the Julian day when 25% of the population arrived were 131-133. Given this data set, it is possible to partially conclude an earlier arrival time as the years continued. Yet, with so many outliers and no clear trend apparent, such a conclusion was not made.
Overall, depending on the species, it seems like the average arrival timings were either around 105-110 or 128-133.
Based on the collection of weather and bird data that underwent data analysis, a linear mixed-effect modeling for each of the species was conducted.
cs.lmer <- lmer(j.day~tmin*tmax*wvec+(1|name),cs.arr.weath,na.action = "na.fail")
Anova(cs.lmer)## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: j.day
## Chisq Df Pr(>Chisq)
## tmin 2.2013 1 0.13790
## tmax 0.2813 1 0.59588
## wvec 0.0249 1 0.87472
## tmin:tmax 1.6650 1 0.19692
## tmin:wvec 6.8876 1 0.00868 **
## tmax:wvec 3.2240 1 0.07257 .
## tmin:tmax:wvec 0.0719 1 0.78864
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1cs.lmer2 <- lmer(j.day~wk.tmin*wk.tmax*wk.wvec+(1|name),cs.arr.weath2,na.action = "na.fail")
Anova(cs.lmer2)## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: j.day
## Chisq Df Pr(>Chisq)
## wk.tmin 0.0330 1 0.855857
## wk.tmax 0.9100 1 0.340111
## wk.wvec 8.3898 1 0.003773 **
## wk.tmin:wk.tmax 8.0046 1 0.004666 **
## wk.tmin:wk.wvec 0.3081 1 0.578849
## wk.tmax:wk.wvec 2.5169 1 0.112631
## wk.tmin:wk.tmax:wk.wvec 0.2901 1 0.590164
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1cs.arr.aic <- dredge(cs.lmer2,fixed = c("wk.tmin","wk.tmax","wk.wvec"),)
cs.kb <- kable(cs.arr.aic[1:4,],caption = "Fit values for nested models of the most complicated lme model")
kable_styling(cs.kb)| (Intercept) | wk.tmax | wk.tmin | wk.wvec | wk.tmax:wk.tmin | wk.tmax:wk.wvec | wk.tmin:wk.wvec | wk.tmax:wk.tmin:wk.wvec | df | logLik | AICc | delta | weight | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 1 | 129.9692 | -0.0110893 | 0.0080211 | 0.0676962 | NA | NA | NA | NA | 6 | -95.70426 | 205.6191 | 0.00000 | 0.9968887699 |
| 2 | 118.3802 | 0.0306095 | 0.1638979 | 0.0944008 | -0.0005833 | NA | NA | NA | 7 | -100.46220 | 217.9514 | 12.33238 | 0.0020926844 |
| 5 | 129.9624 | -0.0113421 | 0.0086465 | 0.0730771 | NA | NA | -4.39e-05 | NA | 7 | -101.81584 | 220.6587 | 15.03965 | 0.0005405406 |
| 3 | 130.0398 | -0.0096051 | 0.0035197 | -0.0285697 | NA | 0.0004057 | NA | NA | 7 | -101.93878 | 220.9046 | 15.28555 | 0.0004780052 |
best.lmer <- lmer(j.day~wk.tmin+wk.tmax+wk.wvec+(1|name),cs.arr.weath2,na.action = "na.fail")
Anova(best.lmer)## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: j.day
## Chisq Df Pr(>Chisq)
## wk.tmin 0.1654 1 0.68427
## wk.tmax 0.8798 1 0.34825
## wk.wvec 4.1757 1 0.04101 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1For the Chimney swift species, it was found that tmin and wind velocity are significant. Through the best model ANOVA analysis, wind velocity is a significant variable.
rh.lmer <- lmer(j.day~tmin*tmax*wvec+(1|name),rh.arr.weath,na.action = "na.fail")
Anova(rh.lmer)## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: j.day
## Chisq Df Pr(>Chisq)
## tmin 0.3174 1 0.57320
## tmax 0.3677 1 0.54425
## wvec 0.3166 1 0.57365
## tmin:tmax 1.6240 1 0.20254
## tmin:wvec 4.7990 1 0.02848 *
## tmax:wvec 3.1283 1 0.07694 .
## tmin:tmax:wvec 0.0741 1 0.78547
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1rh.lmer2 <- lmer(j.day~wk.tmin*wk.tmax*wk.wvec+(1|name),rh.arr.weath2,na.action = "na.fail")
Anova(rh.lmer2)## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: j.day
## Chisq Df Pr(>Chisq)
## wk.tmin 0.7360 1 0.3909
## wk.tmax 0.5582 1 0.4550
## wk.wvec 1.1881 1 0.2757
## wk.tmin:wk.tmax 1.1825 1 0.2769
## wk.tmin:wk.wvec 0.9305 1 0.3347
## wk.tmax:wk.wvec 1.0438 1 0.3069
## wk.tmin:wk.tmax:wk.wvec 0.1270 1 0.7216rh.arr.aic <- dredge(rh.lmer2,fixed = c("wk.tmin","wk.tmax","wk.wvec"),)
rh.kb <- kable(rh.arr.aic[1:4,],caption = "Fit values for nested models of the most complicated lme model")
kable_styling(rh.kb)| (Intercept) | wk.tmax | wk.tmin | wk.wvec | wk.tmax:wk.tmin | wk.tmax:wk.wvec | wk.tmin:wk.wvec | wk.tmax:wk.tmin:wk.wvec | df | logLik | AICc | delta | weight | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 1 | 130.1655 | -0.0063707 | 0.0170196 | -0.0401008 | NA | NA | NA | NA | 6 | -99.05517 | 212.3209 | 0.00000 | 0.9983169742 |
| 5 | 129.8700 | -0.0098545 | 0.0278846 | 0.0499096 | NA | NA | -0.0007693 | NA | 7 | -104.80410 | 226.6352 | 14.31436 | 0.0007779379 |
| 2 | 119.7277 | 0.0315061 | 0.1467828 | -0.0324928 | -0.0004815 | NA | NA | NA | 7 | -105.10711 | 227.2413 | 14.92038 | 0.0005745777 |
| 3 | 130.0711 | -0.0070382 | 0.0196661 | -0.0007893 | NA | -0.0001768 | NA | NA | 7 | -105.66011 | 228.3472 | 16.02638 | 0.0003305102 |
best.lmer <- lmer(j.day~wk.tmin+wk.tmax+wk.wvec+(1|name),rh.arr.weath2,na.action = "na.fail")
Anova(best.lmer)## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: j.day
## Chisq Df Pr(>Chisq)
## wk.tmin 0.6324 1 0.4265
## wk.tmax 0.2580 1 0.6115
## wk.wvec 1.9079 1 0.1672For the Ruby-throated Hummingbird species, it was found that tmin and wind velocity are significant. Through the best model ANOVA analysis, wind velocity was the most significant variable.
bk.lmer <- lmer(j.day~tmin*tmax*wvec+(1|name),bk.arr.weath,na.action = "na.fail")
Anova(bk.lmer)## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: j.day
## Chisq Df Pr(>Chisq)
## tmin 0.2984 1 0.58489
## tmax 0.0729 1 0.78711
## wvec 0.9928 1 0.31907
## tmin:tmax 3.7017 1 0.05436 .
## tmin:wvec 1.9312 1 0.16462
## tmax:wvec 0.1771 1 0.67390
## tmin:tmax:wvec 0.7552 1 0.38483
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1bk.lmer2 <- lmer(j.day~wk.tmin*wk.tmax*wk.wvec+(1|name),bk.arr.weath2,na.action = "na.fail")
Anova(bk.lmer2)## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: j.day
## Chisq Df Pr(>Chisq)
## wk.tmin 0.8856 1 0.34668
## wk.tmax 0.5045 1 0.47751
## wk.wvec 5.4136 1 0.01998 *
## wk.tmin:wk.tmax 0.1043 1 0.74676
## wk.tmin:wk.wvec 0.1581 1 0.69087
## wk.tmax:wk.wvec 0.0692 1 0.79247
## wk.tmin:wk.tmax:wk.wvec 0.2406 1 0.62375
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1bk.arr.aic <- dredge(bk.lmer2,fixed = c("wk.tmin","wk.tmax","wk.wvec"),)
bk.kb <- kable(bk.arr.aic[1:4,],caption = "Fit values for nested models of the most complicated lme model")
kable_styling(bk.kb)| (Intercept) | wk.tmax | wk.tmin | wk.wvec | wk.tmax:wk.tmin | wk.tmax:wk.wvec | wk.tmin:wk.wvec | wk.tmax:wk.tmin:wk.wvec | df | logLik | AICc | delta | weight | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 1 | 107.5975 | -0.0104030 | 0.0193873 | -0.0758179 | NA | NA | NA | NA | 6 | -96.45986 | 207.1302 | 0.00000 | 0.9988266088 |
| 5 | 107.6709 | -0.0089520 | 0.0137142 | -0.1136514 | NA | NA | 0.0005234 | NA | 7 | -102.43040 | 221.8878 | 14.75759 | 0.0006236193 |
| 3 | 107.6877 | -0.0102759 | 0.0170643 | -0.1233338 | NA | 0.0002591 | NA | NA | 7 | -102.86318 | 222.7534 | 15.62314 | 0.0004045462 |
| 2 | 108.7599 | -0.0144746 | -0.0110610 | -0.0796779 | 0.0001214 | NA | NA | NA | 7 | -103.88766 | 224.8023 | 17.67210 | 0.0001452256 |
best.lmer <- lmer(j.day~wk.tmin+wk.tmax+wk.wvec+(1|name),bk.arr.weath2,na.action = "na.fail")
Anova(best.lmer)## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: j.day
## Chisq Df Pr(>Chisq)
## wk.tmin 1.3711 1 0.24163
## wk.tmax 0.7675 1 0.38100
## wk.wvec 5.5154 1 0.01885 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1For the Belted-Kingfisher species, it was found that wind velocity is significant. Through the best model ANOVA analysis, wind velocity and tmin are significant variables.
ys.lmer <- lmer(j.day~tmin*tmax*wvec+(1|name),ys.arr.weath,na.action = "na.fail")
Anova(ys.lmer)## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: j.day
## Chisq Df Pr(>Chisq)
## tmin 0.0735 1 0.7864
## tmax 0.0997 1 0.7522
## wvec 0.1291 1 0.7194
## tmin:tmax 0.0023 1 0.9617
## tmin:wvec 0.4248 1 0.5146
## tmax:wvec 0.0358 1 0.8500
## tmin:tmax:wvec 0.7229 1 0.3952ys.lmer2 <- lmer(j.day~wk.tmin*wk.tmax*wk.wvec+(1|name),ys.arr.weath2,na.action = "na.fail")
Anova(ys.lmer2)## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: j.day
## Chisq Df Pr(>Chisq)
## wk.tmin 0.2267 1 0.6340
## wk.tmax 0.0225 1 0.8808
## wk.wvec 0.2849 1 0.5935
## wk.tmin:wk.tmax 0.0098 1 0.9211
## wk.tmin:wk.wvec 0.3268 1 0.5675
## wk.tmax:wk.wvec 0.3615 1 0.5477
## wk.tmin:wk.tmax:wk.wvec 2.5692 1 0.1090ys.arr.aic <- dredge(ys.lmer2,fixed = c("wk.tmin","wk.tmax","wk.wvec"),)
ys.kb <- kable(ys.arr.aic[1:4,],caption = "Fit values for nested models of the most complicated lme model")
kable_styling(ys.kb)| (Intercept) | wk.tmax | wk.tmin | wk.wvec | wk.tmax:wk.tmin | wk.tmax:wk.wvec | wk.tmin:wk.wvec | wk.tmax:wk.tmin:wk.wvec | df | logLik | AICc | delta | weight | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 1 | 109.1124 | -0.0052289 | 0.0157577 | -0.0410682 | NA | NA | NA | NA | 6 | -132.7204 | 279.6514 | 0.00000 | 0.9972259335 |
| 5 | 109.2259 | -0.0042431 | 0.0113135 | -0.0669357 | NA | NA | 0.0003449 | NA | 7 | -137.8333 | 292.6936 | 13.04218 | 0.0014679776 |
| 3 | 108.8084 | -0.0044762 | 0.0186494 | 0.0295369 | NA | -0.0003736 | NA | NA | 7 | -138.2714 | 293.5697 | 13.91836 | 0.0009472407 |
| 2 | 108.8300 | -0.0042849 | 0.0221159 | -0.0404124 | -2.43e-05 | NA | NA | NA | 7 | -139.2420 | 295.5111 | 15.85967 | 0.0003588483 |
best.lmer <- lmer(j.day~wk.tmin+wk.tmax+wk.wvec+(1|name),ys.arr.weath2,na.action = "na.fail")
Anova(best.lmer)## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: j.day
## Chisq Df Pr(>Chisq)
## wk.tmin 0.1824 1 0.6693
## wk.tmax 0.0372 1 0.8470
## wk.wvec 0.3159 1 0.5741For the Yellow-bellied Sapsucker species, it was found that there is no significance in any variable. Through the best model ANOVA analysis, there’s no specific significance for any variable.
For the Common Yellowthroat species, it was found that wind velocity is significant along with both variables tmin and tmax together. Through the best model ANOVA analysis, wind velocity and tmin are significant variables.
Discussion
After analyzing the data we are able to accurately consider the effects to which weather patterns and climate change have affected the trans-gulf migration of many bird species. We specifically took five species into account: Chimney Swift, Ruby-Throated Hummingbird, Yellow-Bellied Sapsucker, Belted Kingfisher, and the common yellowthroat. We first found the estimated arrival times to the breeding ground over a range of years in attempts to determine if they varied year-by-year. The date in which the species arrives has an effect on its food supply and mating success. While we found differing arrival dates for most of the species, they seemed for the most part sporadic from year to year. It did not seem as if there was a definitive trend over sooner or later arrival dates as time went on. We could estimate the average arrival times to fall within two general ranges: either the species arrived around 105-115 days, or they arrived between 128-133 days. This shows that there was some consistent arrival date between species, but not enough to make any sort of correlation. Next we looked at how weather patterns affected the migration of the species. For most species we found that the tmin was significant and did in fact play a role in the migration of the birds, also in a few of the species we found that wind velocity played a significant role, and lastly there was one species in which tmax played a significant role. There was only one species in which the weather played no role at all. This provides similar results to a similar study (Shamoun-Baranes (n.d.)), in which they concluded that clearer weather conditions were proven to relate to more intense and profound migration patterns in comparison to when weather was not as favorable. Our data is significant because it allows us to understand that while climate change may be affecting the behavior of TGMs, it is hard to determine how many years have to go by before the effects become easily visible to scientists. We can conclude that certain variables, as described above, like certain weather conditions did in fact affect the migration roots of bird species coming to MA. Although there are some studies, such as (Robson and Barriocanal (2011)), in which they conclude that climate change has been associated with certain shifts in biological events such as the spring migratory arrival of birds, we have more difficulty concluding that climate change has a direct effect on that as our Julian days were very sporadic from year to year.
Contributions
Steven worked on R-coding, the Methods write-up, and data analysis. Nick worked on the Introduction, Results, and Discussion sections. Tony helped with problem solving during the coding process along with the Discussion section and reference materials. Amanda helped with the coding portion of the project along with the Results and References sections. Amanda and Steven worked on data analysis. All team members carried out editing at the end.